//
//  KHistogram.m
//  Kia
//
//  Created by Andrey on 22/03/2009.
//  Copyright 2009 Karma Software. All rights reserved.
//

#import <QuartzCore/CoreImage.h>
#import "KHistogram.h"
#import "KImage.h"
#import "KPixel.h"
#import "KSimpleUIntList.h"
#import "KSimpleUIntArray.h"

CGFloat GetEuclideanDistanceBetweenPoints(NSPoint p1, NSPoint p2)
{
	CGFloat deltaX = p2.x - p1.x;
	CGFloat deltaY = p2.y - p1.y;
	
	return (CGFloat)sqrt(deltaX * deltaX + deltaY * deltaY);
}

CGFloat GetWeightForHistogramValue(CGFloat value)
{
	return pow(value, 0.25);
}

void AddValueToHistogram(NSUInteger* histogramData, 
						 NSUInteger intervalCount, 
						 CGFloat value, 
						 NSUInteger* maxPixelsPerInterval)
{
	NSUInteger interval = (NSUInteger)(intervalCount * value);
	if (interval >= intervalCount)
		interval = intervalCount - 1;
	
	assert(interval < intervalCount);
	
	histogramData[interval]++;
	
	if (histogramData[interval] > *maxPixelsPerInterval)
		*maxPixelsPerInterval = histogramData[interval];
}

@implementation KHistogram

@synthesize colorContrast;
@synthesize intervalCount;
@synthesize channelHistograms;

- (id) initWithIntervalCount: (NSUInteger)numIntervals andChannelCount: (NSUInteger)numChannels
{
	if (self = [super init])
	{
		intervalCount = numIntervals;
		channelCount = numChannels;
		channelHistograms = 
		(KChannelHistogramArray)malloc(sizeof(KChannelHistogram) * numChannels);
		
		NSUInteger channelDataSize = sizeof(NSUInteger) * intervalCount;
		for (int i = 0; i < numChannels; i++)
		{
			channelHistograms[i].data = (NSUInteger*)malloc(channelDataSize);
			bzero(channelHistograms[i].data, channelDataSize);
			channelHistograms[i].maxPixelsPerInterval = 0;
			channelHistograms[i].mean = 0.0;
			channelHistograms[i].variance = 0.0;
			channelHistograms[i].uniformity = 0.0;
		}
	}
	
	return self;
}

- (void) computeForImage: (KImage*)image
{
	for (int i = 0; i < channelCount; i++)
	{
		bzero(channelHistograms[i].data, sizeof(NSUInteger) * intervalCount);
		channelHistograms[i].maxPixelsPerInterval = 0;
		channelHistograms[i].mean = 0.0;
		channelHistograms[i].variance = 0.0;
		channelHistograms[i].uniformity = 0.0;
	}
	
	ForEachImagePixel(image)
	{
		KRedundantColorVector colorVector = 
		[[image pixelAtX:x y:y] redundantColorVector];
		
		for (int i = 0; i < channelCount; i++)
		{
			AddValueToHistogram(channelHistograms[i].data, 
								intervalCount, 
								colorVector.data[i], 
								&channelHistograms[i].maxPixelsPerInterval);
		}
	}

	square = image.width * image.height;
	
	[self computeStatistics];
}

- (void) computeForSegment: (NSUInteger)segment ofImage:(KImage*)image withCrispSegmentationMap: (KSimpleUIntArray*)crispSegmentationMap
{
	for (int i = 0; i < channelCount; i++)
	{
		bzero(channelHistograms[i].data, sizeof(NSUInteger) * intervalCount);
		channelHistograms[i].maxPixelsPerInterval = 0;
		channelHistograms[i].mean = 0.0;
		channelHistograms[i].variance = 0.0;
		channelHistograms[i].uniformity = 0.0;
	}
	
	ForEachImagePixel(image)
	{
		NSUInteger index = y * image.width + x;
		
		if (crispSegmentationMap.data[index] == segment)
		{
			KRedundantColorVector colorVector = 
			[[image pixelAtX:x y:y] redundantColorVector];
			
			for (int i = 0; i < channelCount; i++)
			{
				AddValueToHistogram(channelHistograms[i].data, 
									intervalCount, 
									colorVector.data[i], 
									&channelHistograms[i].maxPixelsPerInterval);
			}
		}
	}
	
	square = image.width * image.height;
	
	[self computeStatistics];
}

- (void) computeStatistics
{
	// Compute means
	for (int i = 0; i < channelCount; i++)
	{
		for (int j = 0; j < intervalCount; j++)
		{
			CGFloat intervalLowBound = (CGFloat)j / intervalCount;
			channelHistograms[i].mean += intervalLowBound * channelHistograms[i].data[j];
		}
		
		channelHistograms[i].mean /= square;
	}
	
	// Compute variances
	for (int i = 0; i < channelCount; i++)
	{
		for (int j = 0; j < intervalCount; j++)
		{
			CGFloat intervalLowBound = (CGFloat)j / intervalCount;
			channelHistograms[i].variance += 
			(intervalLowBound - channelHistograms[i].mean) * 
			(intervalLowBound - channelHistograms[i].mean) * 
			channelHistograms[i].data[j];
		}
		
		channelHistograms[i].variance /= square;
	}
	
	// Compute uniformity
//	CGFloat pixelCountMean = (CGFloat)square / intervalCount;
	NSUInteger almostEmptyIntervalSpan = 0;
	CGFloat emptinessThreshold = 5e-4; 
	for (int i = 0; i < channelCount; i++)
	{
		almostEmptyIntervalSpan = 0;
		
		for (int j = 0; j < intervalCount; j++)
		{
			CGFloat relativePixelCount = 
			(CGFloat)channelHistograms[i].data[j] / square;
			
			if (relativePixelCount < emptinessThreshold)
				almostEmptyIntervalSpan++;
			else
			{
				if (almostEmptyIntervalSpan > 0)
				{
					channelHistograms[i].uniformity += 
					almostEmptyIntervalSpan * almostEmptyIntervalSpan;
					almostEmptyIntervalSpan = 0;
				}
			}
		}
		
		// Pick up leftovers
		if (almostEmptyIntervalSpan > 0)
		{
			channelHistograms[i].uniformity += 
			almostEmptyIntervalSpan * almostEmptyIntervalSpan;
		}
		
		channelHistograms[i].uniformity /= (CGFloat)intervalCount * intervalCount;
		channelHistograms[i].uniformity = 1.0 - channelHistograms[i].uniformity;
	}

	// Compute color contrast
	KChannelHistogram hueHistogram = channelHistograms[HueImageChannel];
	CGFloat angleStep = 2 * pi / intervalCount;
	CGFloat weightsSum = 0.0;
	CGFloat gravityCenterX = 0.0;
	CGFloat gravityCenterY = 0.0;
	for (int i = 0; i < intervalCount; i++)
	{
		CGFloat angle = i * angleStep;
		CGFloat value = 
		(CGFloat)hueHistogram.data[i] / square;
		
		NSPoint pointOnCircle = NSMakePoint(cos(angle), sin(angle));
		
		CGFloat weight = GetWeightForHistogramValue(value);
		
		gravityCenterX += pointOnCircle.x * weight;
		gravityCenterY += pointOnCircle.y * weight;
		weightsSum += weight;
	}
	
	gravityCenterX /= weightsSum;
	gravityCenterY /= weightsSum;
	
	NSPoint gravityCenter = NSMakePoint(gravityCenterX, gravityCenterY);
	NSPoint geometricCenter = NSMakePoint(0, 0);
	
	colorContrast = 1.0 - GetEuclideanDistanceBetweenPoints(geometricCenter, gravityCenter);
	
	return;
}

- (void) saveToTikzFile: (NSString*)tikzFilename forChannel: (enum KImageChannel)channel
{
	NSString* imageName = 
	[tikzFilename stringByDeletingPathExtension];
	
	NSString* plotDataFilename = 
	[imageName stringByAppendingString:@".table"];
	
	NSString* plotData = [NSString string];
	for (int i = 0; i < intervalCount; i++)
	{
		CGFloat x = i * 3.0 / intervalCount;
		CGFloat y = 
		(CGFloat)channelHistograms[channel].data[i] / 
		channelHistograms[channel].maxPixelsPerInterval;
		
		plotData = 
		[plotData stringByAppendingFormat:@"%.5f\t%.5f\n", x, y];
	}
	[plotData writeToFile:plotDataFilename 
			   atomically:NO 
				 encoding:NSUTF8StringEncoding 
					error:NULL];
	
	NSString* tikzContent = @"\\begin{tikzpicture}[scale = 3, >=stealth]\n";
	// Draw grid
	tikzContent = 
	[tikzContent stringByAppendingString:@"\t\\draw[very thin, color = gray] (-0.1, -0.1) grid (3.1, 1.1);\n"];
	// Draw axis
	tikzContent = 
	[tikzContent stringByAppendingString:@"\t\\draw[->] (-0.1, 0) -- (3.1, 0) node[right] {$c_b$};\n"];
	tikzContent = 
	[tikzContent stringByAppendingString:@"\t\\draw[->] (0, -0.1) -- (0, 1.1) node[above] {$h_b$};\n"];
	tikzContent = 
	[tikzContent stringByAppendingFormat:@"\t\\draw[thick] plot[smooth,ycomb] file {%@};\n", plotDataFilename];
	
	tikzContent = [tikzContent stringByAppendingString:@"\\end{tikzpicture}"];
	[tikzContent writeToFile:tikzFilename 
				  atomically:NO 
					encoding:NSUTF8StringEncoding 
					   error:NULL];
	
}

- (void) saveCircularHueHistogramToTikzFile: (NSString*)tikzFilename
{
	NSString* tikzContent = @"\\begin{tikzpicture}[scale = 2, >=stealth]\n";
	// Draw grid
	tikzContent = 
	[tikzContent stringByAppendingString:@"\t\\draw[step=.5cm,gray,very thin] (-1.9, -1.9) grid (1.9, 1.9);\n"];

	// Draw axis
	tikzContent = 
	[tikzContent stringByAppendingString:@"\t\\draw[->] (-2, 0) -- (2, 0) node[right] {$x$};\n"];
	tikzContent = 
	[tikzContent stringByAppendingString:@"\t\\draw[->] (0, -2) -- (0, 2) node[above] {$y$};\n"];

	// Draw circle
	tikzContent = 
	[tikzContent stringByAppendingString:@"\t\\draw (0,0) circle (1cm);\n"];
	
	// Draw histogram
	KChannelHistogram hueHistogram = channelHistograms[HueImageChannel];
	CGFloat angleStep = 2 * pi / intervalCount;
	CGFloat gravityCenterX = 0.0;
	CGFloat gravityCenterY = 0.0;
	for (int i = 0; i < intervalCount; i++)
	{
		CGFloat angle = i * angleStep;
		CGFloat value = 
		(CGFloat)hueHistogram.data[i] / hueHistogram.maxPixelsPerInterval;
		
		NSPoint pointOnCircle = NSMakePoint(cos(angle), sin(angle));
		
		NSPoint histogramPlotPoint =
		NSMakePoint((1.0 + value) * cos(angle), 
					(1.0 + value) * sin(angle));
		
		tikzContent = 
		[tikzContent stringByAppendingFormat:@"\t\\draw (%.4f, %.4f) -- (%.4f, %.4f);\n",
		 pointOnCircle.x, pointOnCircle.y, histogramPlotPoint.x, histogramPlotPoint.y];
		
		// Drawing hacks mode ON
		value *= hueHistogram.maxPixelsPerInterval;
		value /= square;
		// Drawing hacks mode OFF
		
		gravityCenterX += pointOnCircle.x * value;
		gravityCenterY += pointOnCircle.y * value;
	}
	
	// Draw gravity center
	tikzContent = 
	[tikzContent stringByAppendingFormat:@"\t\\draw (%.4f, %.4f) circle (0.01cm);\n", gravityCenterX, gravityCenterY];
	
	// Draw geometric center
	tikzContent = 
	[tikzContent stringByAppendingString:@"\t\\draw (0,0) circle (0.01cm);\n"];
	
	// Draw distance
	tikzContent = 
	[tikzContent stringByAppendingFormat:@"\t\\draw (0, 0) -- (%.4f, %.4f) node[below] {$\\rho$};\n", gravityCenterX, gravityCenterY];
	
	// Draw clue
	tikzContent = 
	[tikzContent stringByAppendingFormat:@"\t\\draw[->] (0, 0) +(20: 0.9) arc (20: 70: 0.9) node[below] {$c_h$};\n", gravityCenterX, gravityCenterY];
	
	tikzContent = [tikzContent stringByAppendingString:@"\\end{tikzpicture}"];
	[tikzContent writeToFile:tikzFilename 
				  atomically:NO 
					encoding:NSUTF8StringEncoding 
					   error:NULL];
}

- (enum KImageChannel) maxVarianceChannel
{
	enum KImageChannel result = 0;
	CGFloat maxValue = 1e-30;
	
	for (int i = 0; i < channelCount / 2; i++)
	{
		if (channelHistograms[i].variance > maxValue)
		{
			maxValue = channelHistograms[i].variance;
			result = i;
		}
	}
	
	return result;
}

- (void) dealloc
{
	for (int i = 0; i < channelCount; i++)
	{
		free(channelHistograms[i].data);
	}
	
	free(channelHistograms);
	
	[super dealloc];
}

@end
